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Abstract 

Given a function / defined on a bounded domain Q C and a number A'^ > 0, we study the 
properties of the triangulation 7jv that minimizes the distance between / and its interpolation on the 
C , associated finite element space, over all triangulations of at most A'^ elements. The error is studied in 

the norm X = U" for 1 < p < oo and we consider Lagrange finite elements of arbitrary polynomial 
degree m— 1. We establish sharp asymptotic error estimates as N ^ +oo when the optimal anisotropic 
triangulation is used, recovering the results on piecewise linear interpolation [21 El HI], an improving 
the results on higher degree interpolation [9l 1101 111] . These estimates involve invariant polynomials 
applied to the m-th order derivatives of /. In addition, our analysis also provides with practical 
■ strategies for designing meshes such that the interpolation error satisfies the optimal estimate up to a 

' fixed multiplicative constant. We partially extend our results to higher dimensions for finite elements 

, on simplicial partitions of a domain Q C JR''. 
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^ '. 1 Introduction. 

\Q . 1.1 Optimal mesh adaptation 

O 

In finite element approximation, a usual distinction is between uniform and adaptive methods. In the 
latter, the elements defining the mesh may vary strongly in size and shape for a better adaptation to the 
local features of the approximated function /. This naturally raises the objective of characterizing and 
constructing an optimal mesh for a given function /. 

Note that depending on the context, the function / may be fully known to us, either through an 
explicit formula or a discrete sampling, or observed through noisy measurements, or implicitly defined as 
' the solution of a given partial differential equation. 

^ , In this paper, we assume that / is a function defined on a polygonal bounded domain D, C H^. For a 

given conforming triangulation T of f2, and an arbitrary but fixed integer to > 1, we denote by Im,T th^ 
standard interpolation operator on the Lagrange finite elements of degree to — 1 space associated to T. 
Given a norm X of interest and a number > 0, the objective of finding the optimal mesh for / can be 
formulated as solving the optimization problem 

min 11/ - /m,r/llx, 
#(r)<jv 

where the minimum is taken over all conforming triangulations of cardinality N. We denote by T/v the 
minimizer of the above problem. 

Our first objective is to establish sharp asymptotic error estimates that precisely describe the behavior 
of 11/ — Im,Tf\\x as iV — >■ +00. Estimates of that type were obtained in "SjlJKTS] in the particular case 
of linear finite elements (m — 1 = 1) and with the error measured in X = L^. They have the form 
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limsupf^V min ||/ - I^.t/WlA < C\\^\ det(d2/)|||^., - = - + 1, (1) 
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which reveals that the convergence rate is govc^rncxl by the quantity det((i^/)|, which depends non- 
hnearly tlie Hessian cP f . This is heavily tied to the fact that we allow triangles with possibly highly 
anisotropic shape. In the present work, the polynomial degree m — 1 is arbitrary and the quantities 
governing the convergence rate will therefore depend nonlinear ly on the m-th order derivative dJ^ f . 

Our second objective is to propose simple and practical ways of designing meshes which behave similar 
to the optimal one, in the sense that they satisfy the sharp error estimate up to a fixed multiplicative 
constant. 



1.2 Main results and layout 

We denote by the space of homogeneous polynomials of degree m 

:= Spanja;*^?;' : k + 1 = m}. 

For any triangle T, we denote by Im,T the local interpolation operator acting from C°(r) onto Pm-i the 
space of polynomials of total degree m — 1. The image of v G C°{T) by this operator is defined by the 
conditions 

for all points 7 e T with barycentric coordinates in the set {0, T^^'' " '-'-}• ^® denote by 

em,T{v)p := - Im,TV\\LP{T) 

the interpolation error measured in the norm L^iT). We also denote by 



em,T{v)p ■■= \\V - Im,TV\\LP = I ^ e„,T(v) 



P 
P 

TeT 



the global interpolation error for a given triangulation T, with the standard modification if p = 00. 

A key ingredient in this paper is a function defined by a shape optimization problem: for any fixed 
1 < p < 00 and for any tt G H^, we define 

^rn.pU) := I mf^em,T(7r)p. (2) 

Here, the infimum is taken over all triangles of area \T\ = 1. Note that from the homogeneity of tt, we 
find that ^ 

inf em,Ti'^)p = Km,p{TT)A^~^p . (3) 

\T\=A 

This optimization problem thus gives the shape of the triangles of a given area which is at best adapted 
to the polynomial tt in the sense of minimizing the interpolation error measured in L^. We refer to Km,p 
as the shape function. We discuss in §2 the main properties of this function. 

Our asymptotic error estimate for the optimal triangulation is given by the following theorem. 

Theorem 1.1 For any polygonal domain Q C M^, and any function f e C™(0), there exists a sequence 
of triangulations {Tn)n>No! with #{Tn) = such that 



lim sup N " em,TN {f)p < 



m! 



1 m ^ 1 



An important feature of this estimate is the "lim sup" asymptotical operator. Recall that the upper limit 
of a sequence {un)n>No is defined by 

limsupMjv ■— fim sup u„, 

iV Af->oo„>jv 
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and is in general strictly smaller than the supremum supj^yj^^ un- It is still an open question to find an 
appropriate upper estimate for sup^ -/V'"^^em.T]v(/)p when optimally adapted anisotropic triangulations 
are used. 

In the estimate (|4]), the m-th derivative c?™/ is identified to an homogeneous polynomial in II„j: 

m! d^xd^y k\ II 

k+l=m " 

In order to illustrate the sharpness of we introduce a slight restriction on sequences of triangulations, 
following an idea in [3]: a sequence {Tn)n>No of triangulations, such that #{Tn) — N, is said to be 
admissible if 

sup diam(T) < CaN'^/^, (5) 

for some Ca > independent of N. The following theorem shows that the estimate ^ cannot be improved 
when we restrict our attention to admissible sequences. It also shows that this class is reasonably large 
in the sense that ([4]) is ensured to hold up to small perturbation. 

Theorem 1.2 Let fl C IR^ be a compact polygonal domain, and f € C"^{fl). Denote | ^ + ^- For 
all admissible sequences of triangulations {Tn)n>Noj one has 

liuiMN'^emTNiDp > 

For all £ > 0, there exists an admissible sequence of triangulations {'T§)n>No: such that 

liuisnp em.,r^{f)p < 

Note that the sequences (T'^)n>No satisfy the admissibility condition ^ with a constant C^(e) which 
may explode as e — 0. The proofs of both theorems are given in §3. These proofs reveal that the 
construction of the optimal triangulation obeys two principles: (i) the triangulation should equidistribute 
the local approximation error em,T{f)p between each triangle and (ii) the aspect ratio of a triangle T 
should be isotropic with respect to a distorted metric induced by the local value of d™/ on T (and therefore 

anisotropic in the sense of the euclidean metric). Roughly speaking, the quantity ||i^m,p (4nf) lU'Cr) 
controls the local interpolation L^-error estimate on a triangle T once this triangle is optimized with 
respect to the local properties of /. This type of estimate differs from those obtained in [2] which hold 
for any T, optimized or not, and involve the partial derivatives of / in a local coordinate system which 
is adapted to the shape of T. 

The proof of the upper estimates in Theorem 11.21 involves the construction of an optimal mesh based 
on a patching strategy similar to [4 . However, inspection of the proof reveals that this construction 
becomes effective only when the number of triangles N becomes very large. Therefore it may not be 
useful in practical applications. 

A more practical approach consists in deriving the above mentioned distorted metric from the exact 
or approximate data of d'"/, using the following procedure. To any n e we associate a symmetric 
positive definite matrix Ht, G S*^. If z G il and d"^f{z) is close to tt, then the triangle T containing z 
should be isotropic in the metric hj^. The global metric is given at each point z by 

h{z) = s{Tr,)h^^, 7r^=d"/(z), 

where sijTz) is a scalar factor which depends on the desired accuracy of the finite element approximation. 
Once this metric has been properly identified, fast algorithms such as in [571 123 H] can be used to design 
a near-optimal mesh based on it. Recently in |201 [6], several algorithms have been rigorously proved to 
terminate and produce good quality meshes. Computing the map 

7reH™^->/^^e52+, (6) 



d'^f 



d"'f 



■ e. 



Li(n') 
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is therefore of key use in applications. This problem is well understood in the case of linear elements 
(m = 2): the matrix h-^ is then defined as the absolute value (in the sense of symmetric matrices) of the 
matrix associated to the quadratic form tt. In contrast, the exact form of this map in the case m > 3 is 
not well understood. 

In this paper, we propose algebraic strategies for computing the map ([B]) for m = 3 which corresponds 
to quadratic elements. These strategies have been implemented in an open-source Mathematica code 
|25j . In a similar manner, we address the algebraic computation of the shape function Km,p{T^) from the 
coefficients of tt e Hm, when m > 3. All these questions are addressed in §4, 5 and 6. 

In §4, we discuss the particular case of linear {m = 2) and quadratic (m = 3) elements. In this case, 
it is possible to obtain explicit formulas for K„i^p{'K) from the coefficients of tt. In the case m ~ 2, this 
formula is of the form 

K2,p{ax^ + 26a;?/ + cy^) = cr\/|62 - ac\, 

where the constant cr only depends on p and the sign of — ac, and we therefore recover the known 
estimate ([IJ from Theorem 11.11 The formula for m — S involves the discriminant of the third degree 
polynomial cfif. Our analysis also leads to an algebraic computation of the map ([6]). We want to 
mention that a different strategy for the the construction of the distorted metric and the derivation of 
error estimate for finite element of arbitrary order was proposed in [9l. In this approach, the distorted 
metric is obtained at a point z g £7 by finding the largest ellipse contained in a level set of the polynomial 
d"^fz- This optimization problem has connections with the one that defines the shape function in ^ as 
we shall explain in §2. The approach proposed in the present work in the case m = 3 has the advantage 
of avoiding the use of numerical optimization, the metric being directly derived from the coefficients of 

In §5, we address the case m > 3. In this case, explicit formulas for j,(7r) seem out of reach. 
However we can introduce explicit functions Km(7r) which are polynomials in the coefficients of tt, and 
are equivalent to A'„i,p(7r), leading therefore to similar asymptotic error estimates up to multiplicative 
constants. At the current stage, we did not obtain a simple solution to the algebraic computation of the 
map ([6]) in the case m > 3. The derivation of K,„ is based on the theory of invariant polynomials due to 
Hilbert. Let us mention that this theory was also recently applied in |22j to image processing tasks such 
as affine invariant edge detection and denoising. 

We finally discuss in §6 the possible extension of our analysis to simplicial elements in higher dimension. 
This extension is not straightforward except in the case of linear elements m — 2. 



2 The shape function 

In this section, we establish several properties of the function Km,p which will be of key use in the sequel. 
We assume that m > 2 is an integer, and p € [l,oo]. We equip the finite dimensional vector space 
with a norm || • || defined as the supremum of the coefficients 

m 

If 7r(x,2/) = y^aix'y™"', then ||7r|| = max \a,\. (7) 

1=0 — 

Our first result shows that the function Km.p vanishes on a set of polynomials which has a simple algebraic 
characterization. 

Proposition 2.1 We denote by s,n [-yJ + 1 smallest integer strictly larger than m/2. The 
vanishing set of Km.p is the set of polynomials which have a generalized root of multiplicity at least Sm'- 

Krn,pi'^) = ^ 7r(2;, y) = (ax + /3y)^'"TT, for some a, (5 ]R and tt G ]Hm-Sm- 



Proof: We denote by Teq a fixed equilateral triangle of unit area, centered at 0. 

We first assume that n{x,y) = {ax + /Sy)"™?!-. Then there exists a rotation i? e O2 and tt e H„ 
such that 

J2 a,x'y'^-''"^-' \ , 
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Therefore denoting by <j)^ the Hnear transform (j)^{x,y) — R (ex, ^) we obtain 

IItto 0,11 = max |a,|e2-'"-"+2'' < e^'-^-'"'\\%\\ ^ as e ^ 0. 

Consequently 

em,0,(Te<,)(7r)p = em,T,j7r o (f)^)p ->0 as e 0. 

Since | det0e| = 1, the triangles 4'e{Tcq) have unit area, and therefore Km^pij^) — 0. 

Conversely, let tt G Hm\{0} be such that K,n,p{'^) = 0. Then there exists a sequence (T„)„>o of 
triangles with unit area such that em.T„(j^)p — ^ 0. We remark that the interpolation error eT(T^)p of 
TT G Hm is invariant by a translation : z i— >■ z + /i of the triangle T. Indeed n — tt o G Pm-i so that 

IItI" — ^m,T7I'||Lp(Tfc(T)) = Ik o Th — /m,T(7I' o Th)\\Lp{T) = Ik — -^m, T7r|| Lp(T) ■ (8) 

Hence we may assume that the barycenter of r„ is 0, and write r„ = (/'n(7oq), for some linear transform 
0„ with det 0„ = 1. Since em,T^^{-)p is a norm on H^, it follows that tt o — >• 0. 
The linear transform (f>n has a singular value decomposition 

(pn = UnoDnO Vn, whcrc [/„ , K G O2, and D„ = ^ "^^^ ^ , < e„ < 1. 

Since the orthogonal group O2 is compact, there is a uniform constant C such that 

|koi/|| <C|k||, vreH™, ye02. 

Therefore 

Ik o [/„ o D„ll = Ik ° c^n ° A. o K o Kr'll < cWtt o ,/,„|| ^ 0. 

Denoting by aj,„ the coefficient of x'y™"* in tt o f7„, we find that ai^nE^'^™ tends to as n — > +00. In 
the case where i < Sm, this implies that aj,„ tends to as n — >■ +cx) 

Moreover, again by compactness of O2, we may assume, up to a subsequence, that C/„ converges to 
some U £ 02- Denoting by the coefficient of x^y™^'^ in tt o J7, we thus find that Ui = ii i < Sm- This 
This implies that tt oU{x,y) — a;*'"7r(x,y) which concludes the proof. o 



Remark 2.1 In the simple case m — 2, we infer from Proposition \2.1\ that K2.p{tt) = if and only if 
TT is of the form 7r(x,j/) = x'^ up to a rotation, and therefore a one- dimensional function. For such a 
function, the optimal triangle T degenerates to a segment in the y direction, i.e. optimal triangles of 
a fixed area tend to be infinitely long in one direction. This situation also holds when m > 2. Indeed, 
we see in the second part in the proof of Proposition \2.1\ that if tt is a non-trivial polynomial such that 
Km,p{'^) = 0, then £„ must tends to as n — > +00. This .shows that Tn = 4'n{T) tends to be infinitely 
flat in the direction Uey with Cy — (0,1). However, Km.p{TT) ~ does not any longer mean that tt is a 
polynomial of one variable. 

Our next result shows that the function Km,p is homogeneous, and obeys an invariance property with 
respect to linear change of variables. 

Proposition 2.2 For all tt e III„i, X e IR and e C{R^), 

Kra,p{\^) = |A|X„,p(7r) (9) 

K^,p{^o4>) = \deX4>r/^K^^p{TT) (10) 

Proof: The homogeneity property (|9]) is a direct consequence of the definitions of Km,p- In order to 
prove the invariance property (1101) we assume in a first part that det cf) ^ and we define T := 



_ 4>(T) 



V|dot0| 



and j^{z) := 7r( ^1 det (/)|z) = | det 0|'"/27r(z). 
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We now remark that the local interpolant Im.T commutes with linear change of variables in the sense 
that, when (j) is an invertible linear transform, 

for all continuous function v and triangle T. Using this commutation formula we obtain 

em,T(7r o = I det0|"^/Pe,„^0(T)(7r)p 

= |det0r/2e„^tWf 

Since the map T h- > T is a bijection of the set of triangles onto itself, leaving the area invariant, we 
obtain the relation pop when is invertible. When det cf) = 0, the polynomial tt o </> can be written 
{ax + (3y)"^ so that K^^piP o i^) = by Proposition o 

The functions Km,p are not necessarily continuous, but the following properties will be sufRcient for 
our purposes. 

Proposition 2.3 The Junction K„i^p is upper semi- continuous in general, and continuous if m ~ 2 or 
m is odd. Moreover the following property holds: 

If T^n ^ and Km,p{'n:n) — >■ then Km^pij^) = 0. (12) 

Proof: The upper semi-continuity property comes from the fact that the infimum of a family of upper 
semi-continuous functions is an upper semi-continuous function. We apply this fact to the functions 
TT I— >■ e-m,T{T^)p indexed by triangles which are obviously continuous. 

For any polynomial tt S H2, tt — ax^ + 2bxy + cy'^ , we define det n — ac—b^ . It will be shown in §4 that 
K2,pi''^) = "■pi/l det tt\ , where dp only depends on the sign of detyr. This clearly implies the continuity of 
K2,p- We next turn to the proof of the continuity of Km,p for odd m. Consider a polynomial tt G Hm. If 
Km,p{'^) — then the upper semi-continuity of Km^p, combined with its non-negativity, implies that it is 
continuous at tt. Otherwise, assume that Km,p{'^) > 0. Consider a sequence 7r„ S converging to tt, 
and a sequence (/)„ of linear transformations satisfying det0„ = 1, and such that 

lim Ca (T.„)(7r„) = liminf iir„,p(7r*) := lim inf Kra,p{-n*). 

n— >+oo ^ ^ ^' Tr*-^7r r— ^0 || tt — 7r|| <r 

If the sequence 0„ admits a converging subsequence (fin^ — 4>, it follows that 




This asserts that K^ p is lower semi continuous at tt, and therefore continuous at tt since we already 
know that Km,p is upper semi-continuous. 

If (pn does not admit any converging subsequence, then we invoke the SVD decomposition 0„ = 
Un o Dn o Vn, where Un,Vn & O2 and £)„ = diag(£„, ^), where < e„ < 1. (Here and below, we use 
the shorthand diag(a, b) to denote the diagonal matrix with entries a and b) The compactness of O2 
implies that Un admits a converging subsequence — > U. In particular 7r„^ o Un^ converges to tt o [/. 
Therefore, denoting by a^.n the coefficient of in 7r„ o [/„, the subsequence ai^„j. converges to the 

coefficient of x'^y"^^^ in tto t/. Observe also that e„ — > 0, otherwise some converging subsequence could 
be extracted from </>„. Since 60^(7;,^) (7r„) ~ eT^^{'^n°4>n), the sequence of polynomials 7r„o(/)„ is uniformly 
bounded, and so is the sequence 7r„ o t/„ o Dn- Therefore the sequences (ai_„e^/~™)„>o are uniformly 
bounded. It follows that Oi — {) when i < ^. Since m is odd, this implies that tt o U{x,y) — x^'^TT{x,y) 
and Proposition 12.11 implies that K„i p{-K) = which contradicts the hypothesis K^a,p{T^) > 0. 

Last, we prove property (I12[) . The assumption Km_p{TTn) — > is equivalent to the existence of a 
sequence T„ — 4>niTeq) with det(/i„ — 1 such that em.Tn{'^n)p 0. Reasoning in a similar way as in the 
proof of Proposition 12. 1[ we first obtain that 7r„ o (/)„ — 0, and we then invoke the SVD decomposition 
of (j)n to build a converging sequence of orthogonal matrices Un ^ U and a sequence < e„ < 1 such 
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that if tti^n is the coefficient of s'y™ * in 7r„ o [/„, we have ai^n^"^ ™ — > 0. When i < it foUows that 
ai,n — >■ and therefore tt o U{x, y) — x''"^n{x, y). The result follows from Proposition 12. II o 

We finally make a connection between the shape function and the approach developed in 9 . For all 
TT e Urn, we denote by as the level set of |7r| for the value 1, 

K^^{{x,y)^U\ Hx,y)\ < 1}. (13) 

We now define 

Kiin) = sup \E\ . (14) 

\Ee£, ECA^ / 

where the supremum is taken over the set £ of all ellipses centered at 0. The optimization problem 
defining Kf^^ is equivalent to 

inf{deti? : H e S+ andV z E IR^ , {Hz, z) > \TT{z)f^"'}, (15) 

where S*^ is the cone of 2 x 2 symmetric definite positive matrices. The minimizing ellipse E* is then 
given by {{Hz,z) < 1}. The optimization problem described in (jl5l) is quadratic in dimension 2, and 
subject to (infinitely many) linear constraints. This apparent simplicity is counterbalanced by the fact 
that it is non convex. In particular, it does not have unique solutions and may also have no solution. 

Proposition 2.4 On Mm, one has the equivalence 

cKi < K,n,p < CKi 
with constant < c < C independent of p. 



Proof: Let Tcq denote an equilateral triangle of unit area, and B its circumscribed disk. It is easy to 
see that the inscribed disc is B/2. 

We first show that Km,p < CK^^. Let tt e Hm, and let £"„ be a sequence of ellipsoids inscribed in 
and such that \En\ tends to sup{|i?| : E <E £ , E d Ajr} as n — +00. We write i?„ = A„(^„(B), where 
(pn is a linear transform such that det (/)„ — 1 and A„ > 0. We define the triangle r„ = (pniTcq) which 
satisfies |r„| = 1. We then have 

Ik o 0„ - (/„^T„7r) o 0„||2,P(T„,) 
= \\t^ ° (f'n - Im,n^{Tr O (t)n)\\LP{n^) 

< ||7ro0„ -/„_j.^^(7ro(/)„)||ioo(T,,), 

where we have used the commutation formula (jlip . 

Remarking that Im.T^^^ is a continuous operator from lEI„i to P,n-i in the sense of any norm since 
these spaces are finite dimensional, we thus obtain 

^m,p(7r) < Cllk O 0„||ioo(7^^^) 

< Cllk O 0n|U-(S) 
= Ci|k||ioo(0^(B)) 

= ClA~'"|k|jioo(£;^) 



< Ci 



\En\^-"''^ 



\B\ 



where we have used the fact that k| < 1 in En C A^. Letting n — > +00, we obtain that -ft"m,p(7r) < 
CK^{tt) with C = Ci|Br/2. 

We next prove that cKf-^ < K^.p- Let r„ be a sequence of triangles of unit area such that em.T„{T^)p 
tends to Km,pi'^) as n — )■ +00. As already remarked in ([8]) the interpolation error is invariant by 
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translation. We may therefore assume that the triangles r„ have their barycenter at the origin. Then 
there exists linear transforms (/)„ with det 0„ = 1, such that T„ = </>n(2eq). We now write 

lkllL~(0„(S/2)) < \\M\l°-{T„) 

— C'2e„i,T„(7r)p, 

where we have used the fact that || • and eT^^(-)i are equivalent norms on ]H„i, and that em,Tocj(')p 

is an increasing function of p since |Tcq| — 1. By homogeneity, it follows that if A„ :— (C2em,T„ (7i')p)^^^™j 
we have 

lkllL°=(A„0„(S/2)) < A"C2e„,T„(7I')p = 1. 

Therefore the ellipse := A„(/)„(S/2) is contained in A^r, so that 

Letting n —J- +00, we obtain that cK^ < K^^p with c = C2 I ^ 1 . * 

Remark 2.2 Since Km,p and Kf^^ are equivalent, they must vanish on the same set, and therefore Propo- 
sition is also valid for K^. It also easy to see that Kf^^ satisfies the homogeneity and invariance 
properties stated for Km,p in (0) o,nd as well as the continuity properties stated in Provosition \2.S[ 

Remark 2.3 The continuity of the functions K^ p and can be established when m is odd or equal 
to 2, as shown by Proposition \2.S[ but seems to fail otherwise. In particular, direct computation shows 
that Kf{x^y^ — ey'^) is independent of e > and strictly smaller than Kf{x^y^). Therefore Kf is upper 
semi- continuous but discontinuous at the point x^y^ G IH4. 

3 Optimal estimates 

This section is devoted to the proofs of our main theorems, starting with the lower estimate of Theorem 
II. 2( and continuing with the upper estimates involved in both Theorem 11.11 and 11.21 

Throughout this section, for the sake of notational simplicity, we fix the parameters m and p and use 
the shorthand 

K = Km,p and ctItt) = em,T{Tr)p- 

For each point z E ft we define 

^""^^ ^ TFT 
to! 

where / S C""(ri) is the function in the statement of the theorems. We denote by 

a;(r) := sup ||7r^ - 7r^'||, 

I 2 — 2' II <r 

the modulus of continuity of z i~> ttz with the norm || • || defined by ([T]). Note that w(r) — as r — ?> 0. 
3.1 Lower estimate 

In this proof we will use an estimate by below of the local interpolation error. 

Proposition 3.1 Assume that 1 < p < 00. There exists a constant C > 0, depending on f and Q, such 
that for all triangle T C and z E T, 

eriff > KP{n,)\T\'^+^ - C(diamr)"f |T|w(diamr). (16) 
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Proof: Denoting by fiz G fm the Taylor development of / at the point z up to degree m, we obtain 

f{z + u)- ii,{z + w) = m [ {n,+tu{u) - n,{u)){l - t)"'~^dt. 

Jt=o 

and therefore 

II/-M.IIl-(t) < Codiam(T)"^(diam(T)), 

where Cq is a fixed constant. By construction tt^ is the homogenous part of fi^ of degree m, and therefore 
IJ'z — T^z Pm-1- It follows that for any triangle T, we have 

— Im,TlJ-z = TTz — Im,Tl^z- (17) 

We therefore obtain 

I Gt (/)- Gt (71-^)1 < \\{f - Im,T.f) ~ [t^z - I-m,TT^z)\\LP(T) 

< \T\'^mf-UTf) - (/i. -/™.TM.)llL=»(r) 

= |r|i/^||(/-VT)(/-A^.)IU~(T) 

< C^\T\'/P\\f ^ ti.h^^T) 

< CoCi|T|i/Pdiam(T)"a;(diam(r)) 

where Ci is the norm of the operator / — Im,T in L°°{T) which is independent of T. 
From ^ we know that ctItTz) > \T\^^p K{'Kz), and therefore 

eT(/) > if(7r^)|r|^ + ? - CoCi|T|i/Pdiam(T)'"w(diam(r)). 

We now remark that for all p E [1, oo) the function r is convex, and therefore if a, b, c are positive 

numbers, and a > b — c then > max{0, b — c}^ > b^ — pcbP~^. Applying this to our last inequality we 
obtain 

eriff > i^P(7r^)|r|'¥+i -pCoCi(-ft:(7r^))P-i|T|(f-i)('^ + p)+p diam(r)"cj(diamr). 

Since jTl^P-^^^^ + p^+p = |T|(p-i)t|T| < (diamr)"(P-i)|T|, this leads to 

erifT > KP{n,)\T\'^+^ - C{diamT)"'P\T\uj{diamT), 

where C := pCoCi{sup^^Q K{Trz))P^^ . o 

We now turn to the proof of the lower estimate in Theorem 11.21 in the case where p < oo. Consider 
a sequence (T/v)jv>a'o of triangulations which is admissible in the sense of equation ([5]). Therefore, there 
exists a constant Ca such that 

diamT < CaN^^^^, N > Nq, T € Tn 
For T e Tn, we combine this estimate with (fT6|) . which gives 

Averaging over T, we obtain 

eTifr> I KP{^z)\T\'¥dz-\T\N-'^CYCiu{CAN-''^). 
Jt 

Summing on all T £ Tn, and denoting by the triangle in Tzv containing the point z G fJ, we obtain 
the estimate 

erAfr> I if(7r.)|Tf|^dz-iV-^e(7V), (18) 
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where e{N) := |17|CX^Cw(CaA^"^/^) -> as TV ^ +00. The function z ^ \T^\ is hnked with the 
number of triangles in the foUowing way: 



dz 

WF\ 



On the other hand, with ^ = ^ + i, we have by Holder's inequality, 



K'i(TTz)dz < 



KP(7T^\T,^\'^dz 



ilv 



1 



dz 



l-q/p 



(19) 



Combining the above, we obtain a lower bound for the integral term in (|18p which is independent of Tn'- 



f KP{TT,)\T^\'^dz> ( [ K'^{n,)dz] iV-"f/^ 



Injecting this lower bound in (dH) we obtain ernUY > {In K''{TTz)dzy^'^ - e{N) 
us to conclude 



liminfiV^er„(/) > 

N—^-\-oo 



K1{TTz)dz 



(20) 



which is the desired estimate. 

The case p = 00 follows the same ideas. Adapting Proposition 13. 1[ one proves that 

erif) > K{tt,)\T\'^ - C(diamT)'"w(diamT). 

and therefore 

er«(/) > ||if(7r,)|Tf l^ll^^^^j -7V-¥£(7V), 
where e{N) := C'^Cuj{CaN^^) — >■ as +00. The Holder inequality now reads 



(21) 



K{iTz)-dz < 



Li(f2) 



equivalently 



|if(^,)|Tf|^l 



> 



K{-Kz)-dz] iV"^ 



Combining this with (|2ip , this leads to the desired estimate ([20]) with p = 00 and q — -^- 

Remark 3.1 This proof reveals the two principles which characterize the optimal triangulations. Indeed, 
the lower estimate \20\) becomes an equality only when both inequalities in H6\) and fil9(l are equality. 
The first condition - equality in hi 6^) - is met when each triangle T has an optimal shape, in the sense 
that eriiTz) = Ar(7r^)|T|^ + p for some z Cz T. The second condition - equality in I119\) - is met when the 
ratio between KP{'iTz)\T^ and \T^\~^ is constant, or equivalently K{nz)\T\~^^~p is independent of the 
triangle T . Combined with the first condition, this means that the error cxifY equidistributed over the 
triangles, up to the perturbation by (diamT)'"P|r|a;(diamT) which becomes neglectible as N grows. 



3.2 Upper estimate 

We first remark that the upper estimate in Theorem 11.21 implies the upper estimate in Theorem 11.11 by 
a sub-sequence extraction argument: if the upper estimate in Theorem 11.21 holds . then for all n > there 
exists a sequence {Ti!})n>No such that 



limsup A-er^.(/) < 



K 



d"'f 



Li 



1 

n 



10 



with i = i + ^. We then take Tn ^Tj^ 



n{N) 



where 



n{N) = max <n< N ; iV^er"(/) < 



K 



ml 



Li 



For N large enough this set is finite and non empty, and therefore n{N) is well defined. Furthermore 
n{N) — > +00 as iV — > +00 and therefore 



limsup(iV^er„(/)) < 



K 



LI 



We are thus left with proving the upper estimate in Theorem 11.21 We begin by fixing a (large) number 
M > 0. We shall take the limit A/ — > cx) in the very last step of our proof. We define 

Tf^i = {T triangle, \T\ = 1, bary(r) = and diam(r) < M}, 

the set of triangles centered at the origin, of unit area and diameter smaller than M . This set is compact 
with respect to the Hausdorff distance. This allows us to define a "tempered" version of if = Km,p that 
we denote by Km'- 

Km{t^) = inf eT(7r). 

Since Tm is compact, the above infimum is attained on a triangle that we denote by Tm{t^)- Note that 
the map tt i— )■ Tm{t^) need not be continuous. It is clear that Km{t^) decreases as M grows. Note also 
that the restriction to triangles T centered at is artificial, since the error is invariant by translation 
as noticed in ([5]). Therefore Km{t^) converges to K{tt) as M ^ +00. Since Tm is compact, the map 
TT I— >■ maxTgTM sriT^) defines a norm on H™, and is therefore bounded by Cm||7''|| for some Cm > 0. One 
easily sees that the functions tt i->- eT(7r) are uniformly CAf-Lipschitz for all T € Tm, and so is Km- 

We now use this new function Km to obtain a local upper error estimate that is closely related to the 
local lower estimate in Proposition 13. II 

Proposition 3.2 For zi e 57, let T be a triangle which is obtained from Tm{t^zi) by resettling and 
translation ( T = tTM{T^zi) + zq). Then for any Z2 G T , 



erif) < [KmM + BMLo{max{\zi - Z2\,diam{T)})]\T\ 



(22) 

where Bm > is a constant which depends on M . 
Proof: For all zi, Z2 € il, we have 

STm (tt,^ ) (7^22 ) < (^Tm (tt.i ) (^zi ) + I - 7r^2 1 1 

— Km (t^zi ) + Cm\\ T^zi — TTz^ 1 1 , 

< Km{tTz2) + ^CmWt^Zi - TT^JI, 

< Km{tTz2) +'^C:Mi^{\zi - Z2\)- 

Therefore, if T is of the form T = tTMiT^zi) + zq, we obtain by a change of variable that 

eT(^.J < (KAiiTTz,) + 2Cm^{\zi - Z2I)) |r|^ + i^ 

Let Hz E Pm be the Taylor polynomial of / at the point z up to degree to. Using Equation ([T7| we 
obtain 

erif) < er (Ai^2 ) + eT (/ - /^22 ) 
= eri-r^z^) + erif - IJ-Z2) 

< [Km{t:z2) + 2CMt^i\zi - Z2I)) \T\^+^ + erif - fiz^) 
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By the same argument as in the proof of Proposition 13. 1[ we derive that 

erif - ^^2) < C\T\'^ diam(r)"a;(dianir), 

and thus 

erif) < (i^Af(7r.J + 2CaM\zi - Z2I)) |r|^ + t + C\T\-p diam(T)™L.(diamr). 
Since T is the scaled version of a triangle in Tm, it obeys diam(T)^ < M^|r|. Therefore 
erif) < {Km{tt.,) + {2Cm + CAf™)L.(max{|zi - z^l, diam(T)})) |r|^ + ?, 
which is the desired inequahty with Bm '■— 2Ci\i + CM"\ O 
For some r > to be specified later, we now choose an arbitrary triangular mesh 7?. of fi satisfying 

r > sup diam(i?). 

Ren 

Our strategy to build a triangulation that satisfies the optimal upper estimate is to use the triangles R 
as macro-elements in the sense that each of them will be tiled by a locally optimal uniform triangulation. 
This strategy was already used in [3]. 
For all R G TZ we consider the triangle 

{KMiiTb^) + 2BM^{r))-iTM{TTba). 

which is a scaled version of TMiT^ba) where hn is the barycenter of R. We use this triangle to build a 
periodic tiling Vr of the plane: there exists a vector c such that Tr U forms a parallelogram of side 
vectors a and h, with = c — Tn. We then define 

'Pr ■= {Tr + ma + nb : m,ne 72?} U {T^ + ma + n6 : m, n G 72?}. (23) 

Observe that for all tt e Hm, and all triangles T, T' such that T' — — T one has eT(7r) = eT'(7r) since tt 
is either an even polynomial when m is an even integer, or an odd polynomial when m is odd. Since we 
already know that eT(7r) is invariant by translation of T, we find that the local error exij^) is constant 
on aU T &Vr. 

We now define as follows a family of triangulations Ts of the domain fi, for s > 0. For every R E TZ, 
we consider the elements T D R for T £ sVr, where sVn^ denotes the triangulation Vr scaled by the 
factor s. Clearly {TOR, T £ sVr, R E TZ} constitute a partition of il. In this partition, we distinguish 
the interior elements 

r/"s :^ {T E sVr : T E mt{R) , i? G 7^}, 

which define pieces of a conforming triangulation, and the boundary elements T n i? for T G sVr such 
that T n dR ^ 0. These last elements might not be triangular, nor conformal with the elements on 
the other side. Note that for s > small enough, each R E TZ contains at least one triangle in 77°^, 
and therefore the boundary elements constitute a layer around the edges of TZ. In order to obtain a 
conforming triangulation, we proceed as follow: for each boundary element Tn R, we consider the points 
on its boundary which are either its vertices or those of a neighboring element. We then build the 
Delaunay triangulation of these points, which is a triangulation oi T O R since it is a convex set. We 
denote by 7^'^'^ the set of all triangles obtained by this procedure, which is illustrated on Figure 1. 
Our conforming triangulation is given by 

As s — ?> 0, clearly 

^iT}"") < a^s-' and 1^1 ^ ^bd5. 
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Figure 1: a. An edge (Thick) of the macro-triangulation TZ separating to uniformly paved regions (Tr 
is thick, Vb. is dashed), b. Additional edges (dashed) are added near the interface in order to preserve 
conformity, c. The sets of triangles 77°^ (gray) and T^'^'* (white) 

for some constant Cbd which depends on the macro-triangulation TZ. We do not need to estimate Cbd 
since TZ is fixed and the contribution due to Cbd in the following estimates is neglectible as s — >■ 0. We 
therefore obtain that the number of triangles in T^'^'* is dominated by the number of triangles in Ts °^. 
More precisely, we have the equivalence 

#(7;) ^ #(rr^) ^ J2 = E l^l(^M(vr,J + 2i?MC.(r))^ (24) 



in the sense that the ratio between the above quantities tends to 1 as s 0. The right hand side in (|24 
can be estimated through an integral: 



^#(CT < ^ |i?|(i^Af(vrbJ + 2i?MC.(r))« 
Ren 

Ren •'^ 

< I (^m(^.) + Cnh, - TTbnW + 2BMu:{r)ydz 
Ren-'^ 

< [ {Km{^z) + {2Bm + CMMrWdz 

JQ. 

Therefore, since Cm < Bm, 

#{Ts) < s-^ ( f (KMiTTz) + 3BM^{r)ydz + Q^^s] . (25) 



n 



Observe that the construction of Ts gives a bound on the diameter of its elements 

sup diam(T) < sCa, Ca ■— maxdiam(rfi') 
TeTs Ren 

Combining this with (|24p . we obtain that 

sup diam(r) < CAi^Ts)"^^^ for all s > 
Te% 

which is analogous to the admissibility condition ([5]). 

We now estimate the global interpolation error ||/ — Im,Tsf\\Lp '■= (J^reT sr(/)^)^, assuming first that 
1 < p < oo. We first estimate the contribution of 7^'^'^, which will eventually be neglectible. Denoting 
I'z G Pm-i the Taylor polynomial of / up to degree m — 1 ed, z we remark that 

\\f-Im.Tf\\L^(T) - \\{I-ImM~^b^)\\L^{T) < Ci || / - i^b, |Uoo (j.) < CqCi diam(r)™. 



13 



where Ci is the norm of / — Im,T in L°°{T) which is independent of T and Co only depends on the L°° 
norm of d™/. Remarking that erif) = ||/- Irn,Tf\\Lp{T) < \T\'i'\\f ~ Im,Tf\\L-^{T), we obtain an upper 
boimd for the contribution of 7^'''^ to the error: 



< C^CfCbdS sup diani(T)™f 

with = C^CfC^PCbd- We next turn to the the contribution of T^''^ to the error. If T G T;''''^, 
T G R € TZ, we consider any point zi = z e T and define Z2 = bn the barycenter of R. With such 
choices, the estimate (1221) reads 



erif) < (i^M(7r,) +BMw(max{r,CAs})) |T|- + p. 

We now assume that s is chosen smaU enough such that Cas < r. Geometrically, this condition ensures 
that the "micro-triangles" constituting Ts actually have a smaller diameter than the "macro-triangles" 
constituting TZ. This implies 

eriff < {Km{^.) + BMLo{r)f |T|^+i (26) 
Given a triangle T £ 77°^, T C R £ TZ, and a point z G T, one has 

|T| = s^KM{7rb^) + 2BMi^ir)r'' 

< {Km{tTz) - Cm\Wz - TTb^ll + 2SMC^(r))"« 

< s^KM{TTz) + {2BM~CM)^{r)r'' ■ 

Observing that Bm > Cm, and that p — = q, wc inject the above inequality in the estimate 
which yields 

erifr < /"^ (KMiTTz) + BMUj{r)f \T\. 

Averaging on z g T, we obtain 

erUY < s'"^' / {KM{nz) + BMUj{r)y dz. 

JT 

Adding up contributions from all triangles in 7^, we find 
Combining this with the estimate f|25p we obtain, 



er.#(r.)"/2 < {J^ {KMiTTz) + BMu:{r)Y dz + Q^s^ ' (^^ (i^M(vr.) + iBMO^{r)Y dz + Cbd^j 
and therefore, since - = ^ + - , 

' q 2 p ' 

lim^sup(#(r.)™/'er.) < (^^^ {Km{tIz) + 3BAfC^(r))« dz^ ' . 
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It is now time to observe that for fixed M, 




hm / {KM{TT,) + 3BAMr)) 




KM{-Kz)dz, 



and that 




hmsup(#(rs)"/2er,) < ( / A'«(7r,)dz)' +£. 



This gives us the announced statement of Theorem 11.21 by defining 



SAT := min{s > : #(r,) < N), 



and by setting Tn — Tsf, ■ 

The adaptation of the above proof in the case p = oo is not straightforward due to the fact that the 
contribution to the error of 7^*"^ is not anymore neglectible with respect to the contribution of T^**^- For 
this reason, one needs to modify the construction of 7^'^'^. Here, we provide a simple construction but for 
which the resulting triangulation 7^ is non-conforming, as we do not know how to produce a satisfying 
conforming triangulation. 

More precisely, we define T^'^'s in a similar way as for p < oo, and add to the construction of T^'"^ a post 
processing step in which each triangle is splitted in 4P similar triangles according to the midpoint rule. 
Here we take for j the smallest integer which is larger than — . With such an additional splitting, 
we thus have 



which remains neglectible compared to s . We therefore obtain 

#(Ts) < s-^ (^j {Km{ttz) + 3BMio{r))^dz + Cbds'/') (27) 

Moreover, if T e 77*^^ and T C R E TZ, we have according to the estimate (|22l) 

erif) < (XM(7rb«) +SA/c^(max{r,CAs})) 1^1^. 
By construction |r| = s^{KM{nb„) + 2BMUj{r))-^/'''. This implies ctI/) < s" when Cas < r. Therefore 

STa(/) — max{e7;r<!E, BT-bd} < s™ max{l, Cj^^js"}. 
Combining this estimate with ((27l) yields 




erbd(/) < CoCi max diam(T)" < Ct.s^ , 



with := CoCiC^. We also have 



3/2 



(#(r.)™/2er. 




2 



lim sup 

s->0 



and we conclude the proof in a similar way as for p < oo. 
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4 The shape function and the optimal metric for Unear and 
quadratic elements 

This section is devoted to linear (m = 2) and quadratic (m = 3) elements, which are the most commonly 
used in practice. In these two cases, we are able to derive an exact expression for Kj^ pin) in terms of 
the coefhcients of tt. Our analysis also gives us access to the distorted metric which characterizes the 
optimal mesh. While the results concerning linear elements have strong similarities with those of [4], 
those concerning quadratic elements are to our knowledge the first of this kind, although |10| analyzes a 
similar setting. 



4.1 Exact expression of the shape function 

In order to give the exact expression of Km,p, we define the determinant of an homogeneous quadratic 
polynomial by 

det(aa;^ + 2bxy + cy^) = ac — b^, 

and the discriminant of an homogeneous cubic polynomial by 

disc(ax^ + bx'^y + cxy^ + dy^) = b^c^ - 4ac^ - 46^d + l^abcd - 21a^d^. 

The functions det on E2 and disc on IH3 are homogeneous in the sense that 

det(A7r) = A^detTr, disc(A7r) = A^discTr. (28) 

Moreover, it is well known that they obey an invariance property with respect to linear changes of 
coordinates 0: 

det(7r o 0) = (det 4>f det tt, disc(7r o ^) = (det 4>f disc tt. (29) 
Our main result relates Km,p to these quantities. 
Theorem 4.1 We have for all tt e IH2, 



^2,p(7r) = cTp (det 7r)vrdet7r|, 

and for all vr S -ff/3 , 

^3,p(7r) = cr* (disc 7r)-y I disc7r|, 
where (Jp{t) and <7*(t) are constants that only depend on the sign oft. 

The proof of Theorem 14.11 relies on the possibility of mapping and arbitrary polynomial tt e III2 such 
that det(7r) 7^ or tt e II3 such that disc(7r) 7^ onto two fixed polynomials 7r_ or 7r+ by a linear change 
of variable and a sign change. 

In the case of III2 , it is well known that we can choose 7r_ = a;^ — y ^ and 7r+ ~ x'^ +y'^. More precisely, 
to all TT e H2, we associate a symmetric matrix Qt^ such that n{z) — ((5^z,z). This matrix can be 
diagonalized according to 

-f/^fn' V' [/e02, Ai,A2€]R. 



A2 

Then, defining the linear transform 

0. :=[/ 

and Att — sign(Ai) e { — 1, 1}, it is readily seen that 

A^TT O 



Ai|-2 

IA2I- 



x^ + y^ if det TT > 
x^ — y^ if det tt < 0. 

In the case of Hs, a similar result holds, as shown by the following lemma. 

Lemma 4.1 Let n e IH^. There exists a linear transform (j)^^ such that 



, x(a;2-3?/2) i/disc7r>0 , . 

""'^-^^ x(x2 + 3y2) z/disc7r<0. (2^) 
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Proof: Let us first assume that tt is not divisible by y so that it can be factorized as 

TT = \{x- riy){x - r2y){x - r^y), 

with A G IR and G C If discvr > 0, then the are real and we may assmne ri < r2 < r^. Then, 
defining 

0^ = A(2disc7r ; / : \ > \ I ■ 

an elementary computation shows that ir o if)^ — x{x'^ — 3y^). If disc tt < 0, then we may assume that ri 
is real, r2 and ra are complex conjugates with Im(r2) > 0. Then, defining 

= A(2disc7r)-V3 f nir2 + - 2r2r, i{r2 - r,)ny^ \ 

an elementary computation shows that tt o 0^ = x{x'^ + 3j/^). Moreover it is easily checked that t/)^ has 
real entries and is therefore a change of variable in R^. 

In the case where tt is divisible by y, there exists a rotation U G O2 such that tt :— n o U is not 
divisible by y. By the invariance property p9p we know that discTr = discTr. Thus, we reach the same 
conclusion with the choice (f).^ := U o cjjj^ . o 

Proof of Theorem I4.H for all tt G 1H2 such that det tt ^ and for all change of variable (j) and A 7^ 0, 
we may combine the properties of the determinant in (|28p and (|29p with those of the shape function 
established in Proposition [2?2] This gives us 

K2,p{^) _ i^2,p(A7r o 0) 



VldetTrl ^|det(A7ro0)| 
Applying this with (p = (p^ and A = A^, we therefore obtain 

'^^ ' ^ ' ' \ K2,p{x^ - y^) if detTT < 0. 

This gives the desired result with ap{t) = K2^p{x'^ + y^) for t > and ap{t) = K2^p{x'^ — y^) for t < 0. In 
the case where detTT = 0, then tt is of the form 7r(x,y) — \{ax + fSy)"^ and we conclude by Proposition 
Othat K2,pi'K) = 0. 

For all TT G H3 such that discTT 7^ 0, a similar reasoning yields 

^ ^ 1/TT- TioQ-i / ^3.p(a;(a;2 - 3y^)) disc7r>0, 

K,^,{.) = x/ldiscTTllOS ^ I ^^■^(^(^2 + 3y2)) if disc. < 0. ' 

where the constant 108 comes from the fact that disc{x(x^ — 3y^)) = — disc(a;(a;^ — 3j/^)) = 108. This gives 
the desired result with a*{t) = 108'i K3^p{x{x'^ - iy"^)) for t > and cr* (t) = IQS'^ K3^p{x{x^ + Sy^)) for 
< < 0. In the case where discTr = 0, then tt is of the form 7r(x, y) = {ax + l3y)'^('jx + 6y) and we conclude 
by Proposition 12.11 that _ftr3_j,(7r) =0. o 



Remark 4.2 We do not know any simple analytical expression for the constants involved in Up and a*, 
but these can be found by numerical optimization. These constants are known for some special values of 
p in the case m — 2, see for example 1^1. 

4.2 Optimal metrics 

Practical mesh generation techniques such as in [201 [71 [221 [17] are based on the data of a Ricmannian 
metric, by which we mean a field h of symmetric definite positive matrices 

a; G i-> h{x) G 5*2^. 
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Typically, the mesh generator takes the metric h as an input and hopefully returns a triangulation Th 
adapted to it in the sense that all triangles are close to equilateral of unit side length with respect to this 
metric. Recently, it has been rigorously proved in j24l [6] that some algorithms produce bidimensional 
meshes obeying these constraints, under certain conditions. This must be contrasted with algorithms 
based on heuristics, such as [2 6) in two dimensions, and |27j in three dimensions, which have been 
available for some time and offer good performance |8] but no theoretical guaranties. 

For a given function / to be approximated, the field of metrics given as input should be such that 
the local errors are equidistributed and the aspect ratios are optimal for the generated triangulation. 
Assuming that the error is measured in X = and that we are using finite elements of degree to — 1, 
we can construct this metric as follows, provided that some estimate of tTz = ^f^-* is available all points 
z G ^l. An ellipse such that {E^l is equal or close to 

sup \E\ (31) 

is computed, where A^r^ is defined as in (jTS]). We denote by hy^^ G S2 the associated symmetric definite 
positive matrix such that 

Ez^{{x,y) : {x,y)'^hr,^{x,y) <l} . 

Let us notice that the supremum in pip might not always be attained or even be finite. This particular 
case is discussed in the end of this section. Denoting by > the desired order of the error on each 
triangle, we then define the metric by rescaling hj^^ according to 

1 p 1 

h{z) — -^h„ where :— v^p+^ \Ez\ . 

With such a rescaling, any triangle T designed by the mesh generator should be comparable to the ellipse 
z -\- azEz centered around z the barycenter of T, in the sense that 

z + cia^Ez ^ z + C2azEz, (32) 

for two fixed constants < 2ci < C2 independent of T (recall that for any ellipse E there always exist a 
triangle T such that E cT (l2E). 

Such a triangulation heuristically fulfills the desired properties of optimal aspect ratio and error 
equidistribution when the level of refinement is sufficiently small. Indeed, we then have 

Ikz - Irn,T'^z\\LP(T), 
\T\'^\\tTz - Im.TT^zWL^iT), 
\T\^hz\\L^iT), 
\azEz\p\\nz\\L'=°{a,E,), 

'>TL-{-— . .2^.. 11 

az "lEzlpW-KzU^iE,), 

where we have used the fact that ttz g Hm- 

Leaving aside these heuristics on error estimation and mesh generation, we focus on the main com- 
putational issue in the design of the metric h(z\ namely the solution to the problem (|31]) : to any given 
TT £ Hm, we want to associate ft,^ G S2 such that the ellipse defined by /ijr has area equal or close to 

SUPi5g£,EcA, \E\. 

When m ^ 2 the computation of the optimal matrix ft,^ can be done by elementary algebraic means. 
In fact, as it will be recalled below, /i^r is simply the absolute value (in the sense of symmetric matrices) 
of the symmetric matrix [tt] associated to the quadratic form tt. These facts are well known and used in 
mesh generation algorithms for Pi elements. 

When m > 3 no such algebraic derivation of /i^ from tt has been proposed up to now and current 
approaches instead consist in numerically solving the optimization problem p5|l. see [9]. Since these 



em,T(/)p 
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Figure 2: Maximal ellipses inscribed in A^, tt = x{x'^ — Sy^) or tt = x{x'^ + 3y^). 



computations have to be done extremely frequently in the mesh adaptation process, a simpler algebraic 
procedure is highly valuable. In this section, we propose a simple and algebraic method in the case m — 3, 
corresponding to quadratic elements. For purposes of comparison the results already known in the case 
m = 2 are recalled. 

Proposition 4.2 1. Let tt e IH2 be such that det(7r) 7^ 0, and consider its associated 2x2 matrix 
which can be written as 

Then, an ellipse of maximal volume inscribed in is defined by the matrix 

2. Let TT e if/3 be such that discTr > 0, and cjiT^^ a matrix satisfying Ii30\) . Define 

/i. = (C'f0-i. (33) 



-2/3 



Then h^^ defines an ellipse of maximal volume inscribed in K^^. Moreover det /i^ — — ^ — (discTr) 



3 

3. Let IT G be such that discTT < 0, and (jj^^ a matrix satisfying I130\) . Define 

Then hj^ defines an ellipse of maximal volume inscribed in A-^. Moreover det h^^ = ^ \ disc7r|^ . 

Proof: Clearly, if the matrix defines an ellipse of maximal volume in the set A^r, then for any linear 
change of coordinates </>, the metric (0^^)'^/i7r0~^ defines an ellipse of maximal volume in the set Attoi/). 
When TT E ¥12, we know that A^tt o (/)^ = x"^ + y'^ when detTr > 0, and x'^ — y^ when detTr < 0, where 
IAttI — 1. When tt G H3, we know from Lemma [4.11 that tt o cj)^ = x{x'^ — 3y^) when discTr > and 
x(x^ + 3y^) when discTr < 0. Hence it only remains to prove that when tt e {x^ + y^, x'^ — y'^, x{x^ — 3y^)}, 
then = Id, which means that the disc of radius 1 is an ellipse of maximal volume inscribed in A^r, 
while when tt — x{x'^ + 3y^) we have ft,^ = 2^/"^ Id. 

The case tt = + is trivial. We next concentrate on the case tt = x{x'^ + 3?;^), the treatment of 
the two other cases being very similar. Let E be an ellipse included in A^, tt = x{x'^ + 2>y^). Analyzing 
the variations of the function 7r(cos 0, sin 6*), it is not hard to see that we can rotate E into another ellipse 
E' , also verifying the inclusion E' C A^r, and which principal axes are {x = 0} and {j/ = 0}. We therefore 
only need to consider ellipses of the form kx^ + hy^ < 1. For a given value of h, we denote by k{h) the 
minimal value of k for which this ellipse is included in A^. Clearly the boundary of the ellipse, defined 



19 



by k{h)x'^ + hip = 1, must be tangent to the curve defined by Tr{x,y) = 1 at some point (x, y). This 
translates into the following system of equations 

TT{x,y) = 1, 

hx^ + ky^ = 1, (34) 

kydxTr{x,y) - hxdxTT{x,y) = 0. 

Eliminating the variables x and y from this system, as well as negative or complex valued solutions, 
we find that k{h) — ^37^ when h e (0,2], and k{h) ~ k{2) = 1 when h > 2. The minimum of 
the determinant hk{h) — ^ { j[ + h^) is attained for h = 2i . Observing that k{2^) = 2^ we obtain 
as announced h,^ — 2^1'^ Id and that the ellipse of largest area included in is the disc of equation 
21/3^-2,2 ^_ y2^ <- as illustrated on Figure 2.b. 

The same reasoning applies to the other cases. For tt = — we obtain k{\i) = j^, h £ (0, 00). In this 
case the determinant hk{h) is independent of ft., and we simply choose h = 1 = k{l). For tt = x{x'^ — 3y^) 
we obtain k{h) = when h e (0, 1] and k{h) — fc(l) = 1 when h > \. The maximal volume is attained 
when ft = 1, corresponding to the unit disc, as illustrated on Figure 2. a. o 

Remark 4.3 When tt e if/3 and discTr > a surprising simplification happens : the matrix p3p has 
entries which are symmetric functions of the roots ri,r2,r3. Using the relation between the roots and the 
coefficients of a polynomial, we find the following expression 

If IT = ax^ + ibx^y + 3ca;z/ + d?/^, then ft^r = 2"^3(disc tt)^ 

This yields a direct expression of the matrix as a function of the coefficients. Unfortunately there is no 
such expression when discTr < 0. 

At first sight. Proposition 14.21 might seem to be a complete solution to the problem of building an 
appropriate metric for mesh generation. However, some difficulties arise at points z e where det tt^ = 
or discTT^ =0. If tt € 1H2\{0} and detTr = 0, then up to a linear change of coordinates, and a change 
of sign, we can assume that tt — x'^. The minimization problem clearly yields the degenerate matrix 
ft.^ = diag(l,0), the 2x2 diagonal matrix with entries 1 and 0. If tt S 1H3\{0} and discTr = 0, then up 
to a linear change of coordinates either tt — x'^ or tt = x^y. In the first case the minimization problem 
gives again h^^ — diag(l,0). In the second case a wilder behavior appears, in the sense that minimizing 
sequences for the problem ([3T|) are of the type ft^ — diag(e~^,£^) with e — > 0. The minimization process 
therefore gives a matrix which is not only degenerate, but also unbounded. 

These degenerate cases appear generically, and constitute a problem for mesh generation since they 
mean that the adapted triangles are not well defined. Current anisotropic mesh generation algorithms 
for linear elements often solve this problem by fixing a small parameter (5 > 0, and working with the 
modified matrix ftjr := h.,^ + 51d which cannot degenerate. However this procedure cannot be extended 
to quadratic elements, since h^^y is both degenerate and unbounded. 

In the theoretical construction of an optimal mesh which was discussed in §3.2, we tackled this problem 
by imposing a bound M > on the diameter of the triangles. This was the purpose of the modified 
shape function Km^t^) and of the triangle Tm^t^) of minimal interpolation error among the triangles of 
diameter smaller than AI. We follow a similar idea here, looking for the ellipse of largest area included 
in A^ with constrained diameter. This provides matrices which are both positive definite and bounded, 
and vary continuously with respect to the data tt G IH3. The constrained problem, depending on a > 0, 
is the following: 

sup{|i;| : E e£, E C and diauiE < 2or^l'^}, (35) 

or equivalently 

inf{detiJ : H e s.t. (Hz,z) > \tt{z)\^/"\ z e and H > aid}. (36) 



/ 2(^2 _ ac) be -ad \ 
[ be -ad 2{c^ -bd) r 
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Figure 3: The set A^r (full) and the ellipses E^^^a (dashed) for various values of a > when tt £ H2. 

We denote by E^^^a and /i^^^ the solutions to ([55]) and ([551) . In the remainder of this section, we show 
that this solution can also be computed by a simple algebraic procedure, avoiding any kind of numerical 
optimization. In the case where tt S II2, it can easily be checked that 

[/^...]=f/"f "'"^'n^''"^ n" I ; V (37) 

\ max{|A2|, a|, / 

as illustrated on Figure 3. 

When TT € IH3, the problem is more technical, and the matrix /i7r,Q takes different forms depending 
on the value of a and the sign of disc tt. In order to describe these different regimes, we introduce three 
real numbers < < < /i^r and a matrix C/jr G O2 which are defined as follow. We first define /Itt by 

^,-'/'■.= unn{\\z\\ : |7r(z)| = 1}, 

the radius of the largest disc D^^ inscribed in A^. For z^r such that |7r(z^)| — 1 and H^Trll = ^i-n , we 
define U.„ as the rotation which maps 0^ to the vector (||zTr||, 0). We then define by 

2a;;i/2 := max{diam(i;) : E e £ ; C E C A^}, 

the diameter of the largest ellipse inscribed in A^ and containing the disc ■ In the case where tt is of 
the form (ax + by)^, this ellipse is infinitely long and we set a^^ — 0. We finally define by 

2/3-1/2 _ diam(i;,), 

where E^^ is the optimal ellipse described in Proposition 14.21 In the case where discTr = 0, the "optimal 
ellipse" is infinitely long and we set = 0. It is readily seen that < /^tt < cstt < Mtt- 

All these quantities can be algebraically computed from the coefficients of tt by solving equations of 
degree at most 4, as well as the other quantities involved in the description of the optimal /i7r,a and E^^ ^ 
in the following result. 

Proposition 4.3 For vr G IH^ and a > 0, the matrix h^^.a o,nd ellipse E^y^a o,i"e described as follows. 

1. If a > /iTT, then h-^^a = aid and Ea,7T is the disc of radius oT^I'^ . 

2. If a-rr < a < fiyy , then 

K^o. = Uj (^^^ (38) 

and Ed is the ellipse of diameter IcT^I'^ which is inscribed in A^ and contains D^^. It is tangent 
to 9A^ at the two points z-^ and —z-^. 

3. If < a < then E^^^a is tangent to 9A^ at four points and has diameter 2a~^/'^ . There are 
at most three such ellipses and ii',r.Q is the one of largest area. The matrix /i^_q, has a form which 
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Figure 4: The set (full), the disc i?7r,/j, = -D^ (full), the ellipse i?,r,Q„ (full), and the ellipses Et^^^ 
(dashed) for various values of a > when tt £ H3 and a G (a^,oo). Left: discTr < 0. Right: discTr > 



Figure 5: The set (full), the ellipse -E^.a^ (full), the ellipse -Ett,/?, = E^^ (full), and the eUipses E^^^a 
(dashed) for various values of a > when tt G H3 and a G (/S^ja^). Left: discTr < 0. Right: discTr > 

depends on the sign o/disc7r. 
(i) If discTT < 0, then 

where (p^^ is the matrix defined in Provosition H.'Il and Aq determined by det(/i7r,Q — aid) — 0. 
(a) If disc TT > 0, then 



3\l 



where tjiiT o,nd Aq, are given as in the case disc tt < and where V is chosen between the three 
rotations by 0, 60 or 120 degrees so to maximize |i?a,ir|- 

(Hi) If disc TT = and > 0, then there exists a linear change of coordinates cj) such that nocf) = x^y 
and we have 

V 27AH 



where Aq is determined by det(/iTr,Q — aid) = 0. 
4- If ot < Ptt, then h^y^a — h^r and E^^^a — Et^ is the solution of the unconstrained problem. 

Proof: See Appendix. <> 

Figure 4 illustrates the ellipses E^^^a, a G (a^,oo) when discTr > (4. a) or discTr < (4.b). Figure 5 
illustrates the ellipses E^^^a, a G (/3^,a^) when discTr > (5. a) or discTr < (5.b). Note that when 
a > a,r, the principal axes of £",^,0 are independent of a since t/^ is a rotation that only depends on tt, 
while these axes generally vary when /^tt < a < a^^-, since the matrix (j)^^ is not a rotation. 



22 



Remark 4.4 For interpolation by cubic or higher degree polynomials (m > 4j, an additional difficulty 
arises that can be summarized as follows: one should be careful not to "overfit" the polynomial tt with 
the matrix Ht^. An approach based on exactly solving the optimization problem h31\) might indeed lead to 
a metric h{z) with unjustified strong variations with respect to z and/or bad conditioning, and jeopardize 
the mesh generation process. As an example, consider the one parameter family of polynomials 

Trt = x^y^ + ty'^ e Mi, te [-1,1]. 

It can be checked that when t > 0, the supremum 5+ = sup^^^ bcA, 1^1 is finite and independent oft, but 
not attained, and that any sequence En C Att^ of ellipses such that lim„_j.oo \En\ = •5'+ becomes infinitely 
elongated in the x direction, as n ^ oo. For t < 0, the supremum S- = sup^g^ ecA^ 1^1 independent 

oft and attained for the optimal ellipse of equation x'^ + \t\^^'^y'^ < 1. This ellipse becomes 

infinitely elongated in the y direction ast ~> 0. This example shows the instability of the optimal matrix ft,^ 
with respect to small perturbations of n. However, for all values of t € [—1, 1], these extremely elongated 
ellipses could be discarded in favor, for example, of the unit disc D — {x^ + 2/^ < 1} which obviously 

satisfies D C Ajr^ and is a near-optimal choice in the sense that 2\D\ = S-^ < S- = |£'|a/2(V2 + 1)- 



5 Polynomial equivalents of the shape function in higher degree 

In degrees to > 4, we could not find analytical expressions of K^.p or K^, and do not expect them to 
exist. However, equivalent quantities with analytical expressions are available, under the same general 
form as in Theorem 14.11 the root of a polynomial in the coefficients of the polynomial tt G H^. This 
result improves on the analysis of [llj . where a similar setting is studied. 

In the following, we say that a function R is a polynomial on Hm if there exists a polynomial P of 
TO + 1 variables such that for all (ao, • • • , am) & H™'*'^, 

^(f2<^^^'y"'''^ :-P(ao,--- ,am), 

and we define degR := degP. 

The object of this section is to prove the following theorem 

Theorem 5.1 For all degree m > 2, there exists a polynomial K,„ on Mm, and a constant Cm > such 
that for all tt e Mm, and all 1 <p < oo 



— '■VK,„(^) < Km.A^) < Cm "'^Kmiir), 

where rm = degK^. 

Since for fixed m all functions Km.p, 1 < P < oo, are equivalent on Hm, there is no need to keep track 
of the exponent p in this section and we use below the notation Km — Km,oo- In this section, please do 
not confuse the functions Km and K^, as well as the polynomials Qd and Qd below, which notations are 
only distinguished by their case. 

Theorem 1 5. II is a generalization of Theorem 14. 1[ and the polynomial involved should be seen as a 
generalization of the determinant on 1H2 , and of the discriminant on II3 . Let us immediately stress that 
the polynomial is not unique. In particular, we shall propose two constructions that lead to different 
K-m with different degree r^. Our first construction is simple and intuitive, but leads to a polynomial 
of degree that grows quickly with m. Our second construction uses the tools of Invariant Theory to 
provide a polynomial of much smaller degree, which might be more useful in practice. 

We first recall that there is a strong connection between the roots of a polynomial in H2 or II3 and 
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its determinant or discriminant. 



y i<i<2 J 

disc i X Yl {x - ny) \ = X^{ri - r2f{r2 - rsfirs - nf. 

\ 1<»<3 / 

We now fix an integer m > 3. Observing that these expressions are a "cychc" product of the squares 
of differences of roots, we define 

S(A, ri, • • • ,rm) := A*(ri - r2)^ ■ ■ ■ {rm-i - rm)^{r,n - n)^ . 

Since m > 3, this quantity is not invariant anymore under reordering of the r^. For any positive integer 
d, we introduce the symmetrized version of the d-powers of the cyclic product 

(5d(A,ri, • • • ,r„i) ^ §{X,r^^, ■ ■ ■ ,r^^)'^, 
where S„i is the set of all permutations of {1, • • • , m}. 

Proposition 5.1 For all d > there exists a homogeneous polynomial of degree Ad on Mm, with 
integer coefficients, and such that 

m 

IfT^ = - riy) then Qd (tt) = Qd(A,ri, • • • ,r™). 

1=1 

In addition, Qc( obeys the invariance property 

Qd(7ro,/,) ^ (det0)2""iQrf(7r). (39) 

Proof: We denote by the elementary symmetric functions in the r^, in such way that 
l[ix - r,y) = x^ - cy,x^-^y + a^x^-^y^ -■■■ + (-l)'"f7,„y'". 

4=1 

A well known theorem of algebra (see e.g. chapter IV. 6 in [H]) asserts that any symmetrical polynomial 
in the r^, can be reformulated as a polynomial in the Ui. Hence for any d there exists a polynomial Qd 
such that 

Qd{^,ri, ■ ■ ■ ,r,n) = Qdic^l, ■ ■ ■ jCni). 

In addition it is known that the total degree of Qd is the partial degree of Qd in the variable ri , in our 
case 4d, and that Qd has integer coefficients since Qd has. 

Given a polynomial tt G H„i not divisible by y, we write it under the two equivalent forms 

m 

TT = aox"" + aix"-i?/ + • • • + a^y"" = A Y[{x - r,y). 

i=l 

clearly — X and a-i = (— l)*f^- It follows that 

/-I /\ ^ \4d/=, I \ 4d/\ I —0,1 (— l)™an 
Qd(X,ri, ■ ■ ■ ,rm) = A Qd((7l, ■ ■ ■ ,(Jm) = > ■ ■ • ; 

V ao ao 

Since deg Qd = 4d, the negative powers of ao due to the denominators are cleared by the factor Oq'' and 
the right hand side is thus a polynomial in the coefficients ao, • • • , Um that we denote by Qc;(7r). 
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We now prove the invariance of with respect to hnear changes of coordinates, this proof is adapted 
from |19j . By continuity of Qd, it suffices to prove this invariance property for pairs (tt, </>) such that (j) is 
an invertible hnear change of coordinates, and neither tt or vr o i^^^ is divisible by y. 

Under this assumption, we observe that if vr = A J^™ j^(a; — r^y) and 4* — ( a ) ' then tt o (p^^ — 
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^ JXiLii^ - ^iU) where 

A = A(det 0)-" TT(7 + 5n) and h = , ■ 

2 — 1 

Observing that 

det 4) , . 



it fohows that 

S(A,fi,--- ,f™) = (det0)-^™§(A,ri,... ,r„). 
The invariance property ([39]) foUows readily. 



We now define = 21cm{degQd : 1 < d < ml} where lcm{ai, • • • ,ak} stands for the lowest common 
multiple of {oi, • • • , a^}, and we consider the following polynomial on H™: 

ml 



dog Qd 

Clearly K,„ has degree r„j and obeys the invariance property K„i(7r o (f) ~ (det 0)~'^Km(7r). 
Lemma 5.2 Let tt G -Hm- //K„j(7r) = i/ien XmlTr) = 0. 



Proof: We assume that i^m(7r) ^ and intend to prove that K„j(7r) 7^ 0. Without loss of generality, 
we may assume that y does not divide tt, since K^ij: o U) = Km{T:) and Km(7r o U) — Km(7r) for 
any rotation U. We thus write tt — ^YYiLi{x — riy), where G (D . Since ^^(Tr) 7^ 0, we know from 
Proposition 12.11 that there is no group of s,„ :— [yj + 1 equal roots r^. 

We now define a permutation a* G such that ^^.•(i) 7^ ro-*(i+i) for 1 < i < to — 1 and ro.*(m) 7^ 
ro-*(i). In the case where m = 2m' is even and m' of the are equal, any permutation a* such that 
''(T'(i) — ^(T*(3) = • • • = fo-*(2m'-i) satisfies this condition. In all other cases let us assume that the n 
are sorted by equality : ii i < j < k and — rk then ^ rj ^ rk- If to = 2m' is even, we set 
a*(2i — 1) = i and (T*{2i) = m' + i, 1 < i < m'. If to = 2m' + 1 is odd we set a*{2i) — i, 1 < i < m' 
and a* (2i — 1) = m' + «, 1 < « < to' + 1. For example, cr* (4 1 5 2 6 3 7) when m — 7 and 
cr* = (1 5 2 6 3 7 4 8) when to = 8. With such a construction, we find that |cr*(i) — cr* (i + 1)| > to' if m is 
odd and \a* (i) — a* {i + 1) \ > m' — 1 if m is even, for all 1 < i < m, where we have set cr* (m + 1) :— a* (1). 
Hence a satisfies the required condition, and therefore §(A,rCT.(i), • • • ,fcr*(m)) 7^ 0. 

It is well known that if fc complex numbers ai,-- - ,afe G (D are such that af + • • • + = 0, 
for all 1 < d < A:, then ai = • • • = = 0. Applying this property to the to! complex numbers 
S(A,ro.(i), • • • ,rCT(m)), cr G S„i, and noticing that the term corresponding to a* is non zero, we see that 
there exists 1 < c? < to! such that Q(j(7i') — Qdi^j^i, ■ ■ ■ ,rm) 7^ 0. Since has real coefficients, the 
numbers Qd{T^) are real. Since the exponent r^/ degQ^ is even it follows that Km(7r) > 0, which con- 
cludes the proof of this lemma. o 

The following proposition, when applied to the function Kcq — "^VKm concludes the proof of Theo- 
rem [^3] 

Proposition 5.3 Let m > 2, and let K^q : IHm — > iR+ be a continuous function obeying the following 
properties 

1. Invariance property ; K^qiir o (j)) = | det (/)|'¥_ftrcq(7r). 
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2. Vanishing property ; jor all tt G IHm, i/ A'oq(7r) = then Km{T^) — 0. 
Then there exists a constant C > such that -^K^q < < CK^q on Mm- 

Proof: We first remark that K^q is homogeneous in a similar way as K,n: if A > 0, then applying 
the invariance property to (p ~ A™ Id yields Xoq(7r o (Amid)) = _ftroq(A7r) and |det(/)|^ = A. Hence 
K^qiXir) = Ai^eq(7r). 

Our next remark is that a converse of the vanishing property holds: if K„i{7r) = 0, then there exists 
a sequence </>„ of linear changes of coordinates, det0„ — 1, such that tt o as n oo. Hence 

Kcq{n) = Kcq{TT o 0„) — Xoq(O). Furthermore, K^qiO) = by homogeneity. Hence Kcq{n) = 0. 

We define the set NF^ ■= {tt g H„j : ^^^(Tr) = 0}. We also define a set C by a property 
"opposite" to the property defining NFm. A polynomial vr e ]H,„ belongs to Am if and only if 

IKII < IK ° 4'\\ for a-ll 4' such that det4> — 1- 
The sets NFm and A„i are closed by construction, and clearly NFm H A„, = {0}. We now define 

^m(7r) lim inf K'm(7r') 

r— >0 IItt'— 7r||<r 

the lower semi-continuous envelope of Km- If Kmi"^) = then there exists a converging sequence 7r„ — )• tt 
such that Km{T^n) — > 0. According to Proposition 12.31 it follows that ifm(7r) = and hence tt g iVF,„. 
Therefore the lower semi continuous function Km and the continuous function K^q are bounded below 
by a positive constant on the compact set {tt G A^, ||7r|| = 1}. Since in addition K^q is continuous and 
Km is upper semi-continuous, we find that the constant 



C = sup max 

irGAm,||7r|| = l 



Km{Tr)' Kcqiir) 
is finite. By homogeneity of Km and K^q, we infer that on Am 



^K,q <Km< Km < CK,q. (40) 



Now, for any tt S Hm, we consider tt of minimal norm in the closure of the set {n o (f) : det </) = 1}. 
By construction, we have tt S Am , and there exists a sequence 0„ , det ^„ = 1 such that tt o </)„—>• tt as 
n — >■ cxD. If 71 = 0, then iirm(7r) = i^oq(7r) = 0. Otherwise, we observe that 

Kmiif) < Km{TT) < Kmi-iT) and Kcqiir) = Kcqin). 

Where we used the fact that Km, Km and K^q are respectively lower semi continuous, upper semi con- 
tinuous, and continuous on Km- Combining this with inequality ()40p concludes the proof. o 



A natural question is to find the polynomial of smallest degree satisfying Theorem 15.11 This leads us 
to the theory of invariant polynomials introduced by Hilbert 19 (we also refer to '16 for a survey on 
this subject). A polynomial R on is said to be invariant if /i = ™ "^"g ^ jg a positive integer and for 
all TT € TRm and linear change of coordinates 0, one has 

i?(7ro0) = (det(/))^i?(7r). (41) 

We have seen for instance that K„j and are "invariant polynomials" on 1H,„. 

Nearly all the literature on invariant polynomials is concerned with the case of complex coefficients, 
both for the polynomials and the changes of variables. It is known in particular ^16) that for all m > 3, 
there exists m — 2 invariant polynomials Ri, - ■ ■ Rm-2 on lH„j, such that for any tt (complex coefficients 
are allowed) and any other invariant polynomial R on Hm, 

If i?i(7r) = • • • = i?„i_2(7r) = 0, then R{n) = 0. (42) 
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A list of such polynomials with minimal degree is known explicitly at least when m < 8. Defining 



21cm(degi?i) and K^q := 



, we see that Xoq(7r) = implies K„i(7r) — and hence 



KmiT^) = 0. According to proposition 15. 3[ we have constructed a new, possibly simpler, equivalent of 
K 

For example when m = 2 the list is reduced to the polynomial det, and for m = 3 to the 



polynomial disc. For m 
polynomials 



4, given tt 



ax 



ibx y + 6cx y 



ae 



4bd+3c\ J 



Adxy^ 

b c 
c d 
d e 



ey , the list consists of the two 



therefore Ki{Tx) is equivalent to the quantity Y^|/(7r)|^ + J(7r)2. As m increases these polynomials un- 
fortunately become more and more complicated, and their number m — 2 obviously increases. According 
to [16], for m = 5 the list consists of three polynomials of degrees 4, 8, 12, while for m = 6 it consists of 
4 polynomials of degrees 2, 4, 6, 10. 



6 Extension to higher dimension 

The function K^.p can be generalized to higher dimension d > 2 in the following way. We denote by 
Hm.d the set of homogeneous polynomials of degree m in d variables. For all d-dimensional simplex T, 
we define the interpolation operator Im,T acting from C^iT) onto the space Wm-i,d of polynomials of 
total degree m — 1 in d variables. This operator is defined by the conditions ImTv{"i) = ^(7) for all 
point 7 S r with barycentric coordinates in the set {0, -p;^^, ' ' ' 1 !}■ Following Section §1.2, and 
generalizing Definition we define the local interpolation error on a simplex, the global interpolation 
error on a mesh, as well as the shape function. 
For ah TT € MmA, 

Km.p.di'^) ■= inf Ik - /.,„,T7r||p. 

|T| = 1 

where the infimum is taken on all d-dimensional simplexes T of volume 1. The variant introduced 
in (jl4p also generalizes in higher dimension, and was introduced by Weiming Cao in pi . Denoting by £d 
the set of d-dimensional ellipsoids, we define 

KiAA = ( sup \E\ 

with Att = {2 e M'' : |7r(z)| < 1}. Similarly to Proposition 12.41 it is not hard to show that the functions 
Km,p,d{T^) and j^{Tr) are equivalent: there exists constants < c < C depending only on m,d, such 
that 

'^^mM — Km^p^d < C^m,d- 

Let {Tn)n>o be a sequence of simplicial meshes (triangles if d = 2, tetrahedrons if d = 3, . . . ) of a 
d-dimensional, polygonal open set fi. Generalizing JS]), we say that (7^)„>o is admissible if there exists 
a constant Ca verifying 

sup diam(T) < CaN-^/"^. 

The lower estimate in Theorem 11.21 can be generalized, with straightforward adaptations in the proof. If 
/ G C""(ri) and {Tn)n>No is an admissible sequence of triangulations, then 

liminf Af^em,rN(/)p > 

Where i:=^ + l. 

The upper estimate in Theorem 11.21 however does not generalize. The reason is that we used in its 
proof a tiling of the plane consisting of translates of a single triangle and of its symmetric with respect to 




riijd.p 



LKfl) 
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the origin. This construction is not possible anymore in higher dimension, for example it is well known 
that one cannot tile the space H^, with equilateral tetrahedra. 

The generalization of the second part of Theorem ()1.2p is therefore the following. For all m and d, 
there exists a constant C = C(to, d) > 0, such that for any polygonal open set 51 C H*^ and / e C""(f2) 
the following holds: for all £ > 0, there exists an admissible sequence Tn of triangulations of f2 such that 

limsup7V^e„i,r«(/)p < C 

TV— >-oo 

The "tightness" Theorem 11.21 is partially lost due to the constant C. This upper bound is not new, and 
can be found in [5] . In the proof of the bidimensional theorem we define by a tiling Vr of the plane 
made of a triangle Tr, and some of its translates and of their symmetry with respect to the origin. In 
dimension d, the tiling Vr cannot be constructed by the same procedure. The idea of the proof is to 
first consider a fixed tiling Vq of the space, constituted of simplices bounded diameter, and of volume 
bounded below by a positive constant, as well as a reference equilateral simplex Toq of volume 1. We 
then set Vr — (t>{Vo), where (/) is a linear change of coordinates such that Tr = 0(Tcq). This procedure 
can be applied in any dimension, and yields all subsequent estimates "up to a multiplicative constant" , 
which concludes the proof. 

Since this upper bound is not tight anymore, and since the functions Km.p,d are all equivalent to 
^ as p varies (with equivalence constants independent of p) , there is no real need to keep track of the 
exponent p. We therefore denote by Km,d the function Km,oo,d- 

For practical as well as theoretical purposes, it is desirable to have an efhcient way to compute the 
shape function Km,d, and an efhcient algorithm to produce adapted triangulations. The case m = 2, 
which corresponds to piecewise linear elements, has been extensively studied see for instance pi 112). In 
that case there exists constants < c < C, depending only on d, such that for all tt e W.2,d, 

c^|det7r| < i^2,d(7r) < C^|det7r|. 

where det tt denotes the determinant of the symmetric matrix associated to tt. Furthermore, similarly 
to Proposition 14. 2[ the optimal metric for mesh refinement is given by the absolute value of the matrix 
of second derivatives, see [4j 112) . which is constructed in a similar way as in dimension d = 2: with 
U and D = diag(Ai,--- ,Xd) the orthogonal and diagonal matrices such that [tt] = U^DU and with 
\D\ := diag(|Ai|, • • • , |Ac;|), we set hj^ = U'^\D\U. It can be shown that the matrix /i^ defines an ellipsoid 
of maximal volume included into the set A^. The case m = 2 can therefore be regarded as solved. 

For values (m, d) both larger than 2, the question of computing the shape function as well as the 
optimal metric is much more difficult, but we have partial answers, in particular for quadratic elements 
in dimension 3. Following §5, we need fundamental results from the theory of invariant polynomials, 
developed in particular by Hilbert |19]. In order to apply these results to our particular setting, we need 
to introduce a compatibility condition between the degree m and the dimension d. 

Definition 6.1 We call the pair of numbers m > 2 and d > 2 "compatible" if and only if the following 
holds. For all n G IHm,d such that there exists a sequence ((f>n)n>o of d x d matrices with complex 
coefficients, verifying det0„ = 1 and lim„_j.oo T^°4'n — 0, there also exists a sequence ipn ofdxd matrices 
with real coefficients, verifying det = 1 and lim„_j.oo tt o = 0. 

Following Hilbert [19], we say that a polynomial Q of degree r defined on ¥im,d is invariant if /i = ^ 
is a positive integer and if for all tt S and all linear changes of coordinates 4>, 

g(7ro0) = (det(/.f QW. (43) 

This is a generalization of (|4ip . We denote by Im,d the set of invariant polynomials on Km.d- It is easy 
to see that if tt e ltim,d is such that Km,d{T^) = 0, then Q{tt) ~ for all Q G Im.d- Indeed, as seen in 
the proof of Proposition 12.11 if Krn,d{T^) = then there exists a sequence </>„ such that det(/)„ — 1 and 
TT o 0. Therefore (|43|) implies that Q{tt) = 0. The following lemma shows that the compatibility 

condition for the pair (m, d) is equivalent to a converse of this property. 

Lemma 6.1 The pair (m, d) is compatible if and only if for all vr e IHm.d 

Km,dM ^0 if and only if Q{t:) = for all Q G Em,d- 
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Proof: We first assume that the pair (m, d) is not compatible. Then there exists a polynomial 
TTo G Hm^d such that there exists a sequence (pn, det(/)„ = 1 of matrices with complex coefficients such 
that TT o 0„ — > 0, but there exists no such sequence with real coefficients. This last property indicates 
that Km,d{'^) > 0. On the contrary let Q e Im.rf be an invariant polynomial, and set /i = HL^^Q. The 
identity 

QKo0) = (det0)''Q(7ro) 

is valid for all (f) with real coefficients, and is a polynomial identity in the coefficients of (j). Therefore it 
remains valid if (j) has complex coefficients. If follows that Q{no) — Q{tto o (/)„) for all n, and therefore 
Q{tto) = 0, which concludes the proof in the case where the pair (m,d) is not compatible. 

We now consider a compatible pair (m, d). Following Hilbert [TO], we say that a polynomial tt e H^.d 
is a null form if and only if there exists a sequence of matrices with complex coefficients such that 
det (f>n = 1 and tt o 0„ — > 0. We denote by NFm,d the set of such polynomials. Since the pair {m,d) is 
compatible, note that tt G NFm,d if and only if there exists a sequence (/)„ of matrices with real coefficients 
such that det = ! and tt o 0„ — > 0. Hence, we find that 

(p 

Denoting by ^ the set of invariant polynomials on Hm^^ with complex coefficients, a difficult theorem 
of [in] states that 

iVF„M = e H„^d : Qin) = for aU Q e ij^ ,^1 

It is not difficult to check that if Q = Qi + iQ2 where Qi and Q2 have real coefficients then holds 
for Q if and only if it holds for both Qi and Q2, i.e. Qi and (52 are also invariant polynomials. Hence 
denoting by Im.d the set of invariant polynomials on H^.d with real coefficients, we have obtained that 

NF^,d = {tt e 'Hm,d ■■ Q(^) = for aU Q e I™,4 
which concludes the proof. o 



Theorem 6.1 // the pair {m,d) is compatible, then there exists a polynomial K on IHm,d (we set r = 
degK ) and a constant C > such that for all tt G lHm,d 

^ v/kR < KrnA^) < C v/kR • (44) 
If the pair (m,d) is not compatible, then there does not exist such a polynomial K. 

Proof: The proof of the non-existence property when the pair (m, d) is not compatible is reported in 
the appendix. Assume that the pair (m, d) is compatible. We follow a reasoning very similar to §5 to 
prove the equivalence (|44p . 

We use the notations of Lemma 16.11 and consider the set 

= e H„,,d : K„,,d(vr) = 0} - e K„.,d : Q{tt) = 0, Q G I,„,4. 

The ring of polynomials on a field is known to be Noetherian. This implies that there exists a finite 
family Qi, - ■ ■ ,Qs G Im.d of invariant polynomials on H^.d such that any invariant polynomial is of the 
form '}2,PiQi where Pi are polynomials on 11^, ^. We therefore obtain 

NF^^d = e H™,d : QiW = • • • = Q,(7r) - 0}. 

which is a generalization of ()42|) . however with no clear bound on s. 

We now fix such a set of polynomials, set r := 21cmi<i<s degQi, and define 

s 

i=l 
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Clearly K is an invariant polynomial on H„i^d, and NFm,d = {tt G Km,d ■ K(7r) = 0}. Hence the 
function Kcq is continuous on ltlm,d, obeys the invariance property A'oq(7r o (jj) — \ det (f>\Kcq(Tr), and for 
all TT e Em, A'cq(7r) — implies K(7r) — and therefore Km,d{'^) = 0. We recognize here the hypotheses 
of Proposition 15. 3[ except that the dimension d has changed. Inspection of the proof of Proposition 15.31 
shows that we use only once the fact that d = 2, when we refer to Proposition 12.31 and state that if 
(7r„) e Him, T^n ^ and Km{T^n) — ^ 0, then Km{T^) = 0. This property also applies to Km,d, when the 
pair {m^d) is compatible. Assume that (7r„) e Hm.d, 7r„ — tt and that Km,di'^n) 0. Then there exists 
a sequence of linear changes of coordinates 0„, det 0„ = 1, such that 7r„ o 0„ — >• 0. Therefore 

K(7r) — lim K(7r„) — lim K(7r„ o 0„) = 

It follows that TT S NFm,d, and therefore K„i,d{'^) — 0. Since the rest of the proof of Proposition l5. 31 never 
uses that d — 2, this concludes the proof of Equivalence p4)) . o 

Hence there exists a "simple" equivalent of Km,d for all compatible pairs {m,d), while equivalents 
of Km,d for incompatible pairs need to be more sophisticated, or at least different from the root of a 
polynomial. This theorem leaves open several questions. The first one is to identify the list of compatible 
pairs {m,d). It is easily shown that the pairs (to, 2), to > 2, and {2,d), d>2 are compatible, but this 
does not provide any new results since we already derived equivalents of the shape function in these cases. 
More interestingly, we show in the next corollary that the pair (3, 3) is compatible, which corresponds to 
approximation by quadratic elements in dimension 3. There exists two generators S and T of /s.a, which 
expressions are given in [23] and which have respectively degree 4 and 6. 

Corollary 6.1 y^\S\^ + is equivalent to Ks^s on IHa^s. 

Proof: The invariants S and T obey the invariance properties S'(7ro (p) — (det (P)'^S{tt) and T{tto (p) — 
(det(?!))^T(7r). We intend to show that if tt e Hg^g and S{7r) = T(7r) then -ft:3,3(7r) = 0. Let us first 
admit this property and see how to conclude the proof of this corollary. According to Lemma 16.11 the 
pair (3,3) is compatible. The function K^q '■— \/|S'P + is continuous on Ha^a, obeys the invariance 
property Kcq{TT o (f>) = \ det (t)\Kcq{TT) and is such that Kcqiir) = implies K3^q{tt) = 0. We have seen in 
the proof of Theorem 16.11 that these properties imply the desired equivalence of K^q and ^^3.3. 

We now show that S{tt) — T{ti) — implies Kj,,^{tt) = 0. A polynomial tt e II3.3 can be of two 
types. Either it is reducible, meaning that there exists tti S Hi, 3 (linear) and 1^2 S H2,3 (quadratic) such 
that TT — 7ri7r2, or it is irreducible. In the latter case according to |18| . there exists a linear change of 
coordinates (p and two reals a, b such that 

TT o (j) ~ y^z — {x'^ + 3axz^ + bz^). 

A direct computation from the expressions given in |;23^ shows that S{tt o (p) — a and T{tt o 0) = —Ah. 
If 5(7r) = T{tt) = then ^(Tr o 0) = T(7r o 0) = and tt o = y^z - x'^ . Therefore for ah A 7^ 0, 
TT o (/)(Aa;, A^y, A~^z) = Xy^z — X^x^, which tends to as A — 0. We easily construct from this point a 
sequence </>„, det 0„ = 1, such that tt o </>„ — > 0. Therefore 7^3.3 (vr) ~ 0. 

If TT is reducible, then tt ~ ttitt2 where tti is linear and tt2 is quadratic. Choosing a linear change of 
coordinates such that tti o cj) = z we obtain 

TT o (j) = 3z(ax^ + 2bxy + cy^) + z^(ux + vy + wz), 

for some constants a,b,c,u,v,w. Again, a direct computation from the expressions given in |23| shows 
that ^(Troc/)) = -(ac-62)2 (andT(7ro0) = 8{ac-b'^)^). Therefore if S'(7r) =T{tt) = then the quadratic 
function ax^ + 2bxy + cy^ of the pair of variables (x, y) is degenerate. Hence there exists a linear change 
of coordinates ip, altering only the variables x,y, and reals ii,u',v' such that 

TT o (j) o ^ = ^zx^ + z'^{u'x + v'y + wz). 

It follows that TT o (j) o ipix, X^^y, Xz) tends to as A 0. Again, this implies that K3^3{tt) = 0, and 
concludes the proof of this proposition. o 
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We could not find any example of incompatible pair (m, d), which leads us to formulate the conjecture 
that all pairs {m^d) are compatible (hence providing "simple" equivalents of Km,d in full generality). 
Another even more difficult problem is to derive a polynomial K of minimal degree for all couples (to, d) 
which are compatible and of interest. 

Last but not least, efficient algorithms are needed to compute metrics, from which effective triangu- 
lations are built that yield the optimal estimates. A possibility is to follow the approach proposed in [9], 
i.e. solve numerically the optimization problem 

inf{deti? : H eSj and Vz G R^, {Hz,z) > |7r(z)|2/'"}, 

which amounts to minimizing a degree d polynomial under an infinite set of linear constraints. When 
d > 2, this minimization problem is not quadratic which makes it rather delicate. Furthermore, numerical 
instabilities similar to those described in Remark l44l can be expected to appear. 

7 Conclusion and Perspectives 

In this paper, we have introduced asymptotic estimates for the finite element interpolation error measured 
in LP when the mesh is optimally adapted to the interpolated function. These estimates are asymptotically 
sharp for functions of two variables, see Theorem 11.21 and precise up to a fixed multiplicative constant 
in higher dimension, as described in §6. They involve a shape function Km,p (or Km,d,p if d > 2) which 
generalizes the determinant which appears in estimates for piecewise linear interpolation [l2 | l4lfT3 ] . This 
function can be explicitly computed in several cases, as shows Theorem 14.11 and has equivalents of a 
simple form in a number of other cases, see Theorems 15.11 and 16.11 

All our results are stated and proved for sufficiently smooth functions. One of our future objectives is to 
extend these results to larger classes of functions, and in particular to functions exhibiting discontinuities 
along curves. This means that we need to give a proper meaning to the nonlinear quantity K„i p ^477!^) 
for non-smooth functions. 

This paper also features a constructive algorithm (similar to [1]), that produces triangulations obeying 
our sharp estimates, and is described in §3.2. However, this algorithm becomes asymptotically effective 
only for a highly refined triangulation. A more practical way to produce quasi-optimal triangulations is 
to adapt them to a metric, see [6j[20l[7]. This approach is discussed in §4.2. This raises the question of 
generating the appropriate metric from the (approximate) knowledge of the derivatives of the function to 
be interpolated. We addressed this question in the particular case of piecewise quadratic approximation 
in two dimensions in Theorems 14.21 and 14.31 

We plan to integrate this result in the PDE solver FreeFem-|--|- in a near future. Note that a Mathe- 
matica source code is already available on the web [25] . We also would like to derive appropriate metrics 
for other settings of degree to and dimension d, although, as we pointed it in Proposition 14. 4[ this might 
be a rather delicate matter. 

We finally remark that in many applications, one seeks for error estimates in the Sobolev norms W^'^ 
(or VF™'P) rather than in the norms. Finding the optimal triangulation for such norms requires a new 
error analysis. For instance, in the survey |24| on piecewise linear approximation, it is observed that the 
metric /i^ = (evoked in Equation ([37)) ') should be replaced with h-^ — {d'^f)^ for best adaptation in 

norm. In other words, the principal axes of the positive definite matrix h^^ remain the same, but its 
conditioning is squared. 

APPENDIX 

A Proof of Proposition 14.3! 

We consider a fixed polynomial tt e H3, a parameter a > 0, and look for an ellipse E.^ ^ of maximal 
volume included in the set a^^/^UH A^. Since this set is compact, a standard argument shows that there 
exists at least one such ellipse. 

If a > ^TT, then a'^^^D C and therefore a~^/'^D n A^r = a~^^^D. It follows that E^^^a = a'^^'^D, 
which proves part 1. 
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In the following we denote by E'^ the ellipse defined by the matrix (p8|) . Note that any ellipse 
containing Dt, and included in A^r must be tangent to (?A^ at the point Zjr, and hence of the form E'^ for 
some ^ > 0. Clearly Eg C E'^ if and only ii 6 > ii. Therefore E'^ C A^r if and only ii a > a-^- Let E be an 
arbitrary ellipse, let Di the largest disc contained in E, and D2 the smallest disc containing E. Then it 
is not hard to check that \E\ — y/\Di\\D2\- For any a verifying ctTv < a < ^t^, the ellipse E'^ is such that 
Di = -Dtt, which is the largest centered disc contained in A^r, and D2 — a^^^^D, which corresponds to 
the bound 2q!~^/^ on the diameter of Eyr^a- It follows that E'^ is an ellipse of maximal volume included 
in n Ajr, and this concludes the proof of part 2. 

Part 4 is trivial, hence we concentrate on part 3 and assume that (3^^ < a < a^r. 

An elementary observation is that iJ^.a must be "blocked with respect to rotations" . Indeed assume 
for contradiction that Re{ET^^a) C A^ for S [0, e] or [— e, 0], where we denote by Rq the rotation of angle 
6. Observing that the set U0g[o^£]i?6i(-E'7r,Q) contains an ellipse of larger area than E^ ^ and of the same 
diameter, we obtain a contradiction. 

In the following, we say that an ellipse E is quadri-tangent to A^r, when there are at least four points 
of tangency between dE and SAjr (a tangency point being counted twice if the radii of curvature of dE 
and cJAtt coincide at this point). 

The fact that iS^r.a is "blocked with respect to rotations" implies that it is either quadri-tangent to 
Att or tangent to dh-,^ at the extremities of its small axis. In the latter case the extremities of the small 
axis must clearly be the points z-n and — z^r, the closest points of 9A,r to the origin. It follows that -Ett.q 
belongs to the family E'^, 5 > described above, and therefore is equal to E'^ since a < a.^. But E'^ 
is quadri-tangent to A^r, since otherwise we would have E'^ _^ C A^ for some £ > 0. 

We have now established that -E^.q is quadri-tangent to Ajr when /?7r < < Q^ir- This property is 
invariant by any linear change of coordinate: if an ellipse E is quadri-tangent to Ajro^, then 4>{E) is 
quadri-tangent to A^. Furthermore if E is defined by a symmetric positive definite matrix H, then 4>{E) 
is defined by {(f>^^)'^ H(j)~^ . This remark leads us to identify the family of ellipses quadri-tangent to 9A^ 
when TT is among the four reference polynomials x{x'^ — 3y^), x{x'^ + 3y^), x^y and x^. In the case of 
x^ there is no quadri-tangent ellipse and we have OfTr = 0, therefore part 3 of the theorem is irrelevant. 
In the three other cases, which respectively correspond to part 3 (i), (ii) and (iii), the quadri-tangent 
ellipses are easily identified using the symmetries of these polynomials and the system of equations (j34p . 

The ellipses quadri-tangent to x{x^ + 3y^) are defined by matrices of the form Hx = diag(A, ^p-), 
where < A < 2. Note that detH\ is decreasing on (0,2^] and increasing on [2^,2]. Given tt with 
discTT < 0, the optimization problem (j36|) . therefore becomes 

mm{detHA : {<f>-T Hx<f>-' > aid}. 

If the constraint is met for A = 2^^^, we obtain E-jj^a = E-^ and therefore a < (3-^. Otherwise, using the 
monotonicity of A 1— > dctH\ on each side of its minimum 23 we see that the matrix Hx — a(f)^(t>T,: must 
be singular. Taking the determinant, we obtain an equation of degree 4 from which A can be computed, 
and this concludes the proof of part 3 (i). 

The ellipses quadri-tangent to x{x — 3y2) are defined by Hxy = diag(A, ^^)V , where < A < 1 
and y is a rotation by 0, 60 or 120 degrees. Since det Hx,v is a decreasing function of A on (0, 1], we can 
apply the same reasoning as above to polynomials tt such that discTr > 0. This concludes the proof of 
part 3 (ii). 

Last, the ellipses quadri-tangent to xy'^ are defined by Hx — diag(A, 2^^), A > 0. The determinant 
is a decreasing function of A, with lower bound as A 00, and the same reasoning applies again hence 
concluding the proof of part 3 (iii). 

B Proof of non existence property in Theorem 16.11 

Let (rn, d) be an incompatible pair. We know from Lemma 16.11 that there exists ttq G ^m.d such that 
Km,d{T^o) > and Q{nQ) — for all invariant polynomial Q e 11m, d- 

We assume for contradiction that a polynomial K satisfies inequalities pij) . Up to replacing K with 
K^'', we can assume that K takes non negative values on Km,d and that /i = ^ is an integer. The rest of 
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this proof consists in showing that K needs to be an invariant polynomial, thus leading to a contradiction 
since we would then have K(7ro) ~ 0. For this purpose we derive from inequalities (|44|) . and from the 
invariance of Km 4 with respect to changes of variables, the inequalities 

C-2''(det 4>T'K{-k) < K(7r o 0) < C2''(det 0)''K(7r), (45) 

where C is the constant appearing in inequalities (|44l) . We regard the function Q(7r, 0) = K(7r o 0) as a 
polynomial on the vector space V — Hm.d x Md, where Md denotes the space oidxd matrices, and observe 
that it vanishes on the hypersurface Vdct ■= {(■"', </>) G 1^ : det (p — 0}. Since (j) 1— > det((/)) is an irreducible 
polynomial, as shown in [5], it follows that Q(7r, 0) = (det (/))Qi(7r, 0) for some polynomial Qi on V. 
Injecting this expression in inequality (|45|) we obtain that Qi(7r, 0) also vanishes on the hypersurface Vdot 
and the argument can be repeated. By induction we eventually obtain a polynomial K on such that 
K(7r o (/)) = (det (I))'^K{tt, (/>). It follows from inequality (gS]) that for all (tt, (p) & V 

C-^''K{Tr) < K{Tr,(f>) < C'^''K{n). 

This implies that K(7r, does not depend on cj). Otherwise, since it is a polynomial, we could find 
TTi e Hm.d and a sequence (j)„ € Md such that |K(7ri,0„)| — >■ 00. Therefore 

K(7ro0) = (det (/>)^K(7r, (/)) = (det <?!))''K(7r, Id) = (det 0)''K(7r). 

This establishes the invariance property of K, in contradiction with our first argument, and concludes 
the proof. 
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